## Evaluation of Wahlomat experiment -----

# load packages
source("packages.r")


## data import ------------------------------

load("../data/yougov_panel_ger_wide_wahlomat_prep.RData")



### Model main effects -----------------------------

# setup instrument, treatment, covariates
assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "age.1", "educ.1", "hhincome.1", "polinterest.1", "leftright.1")


## Generic function to identify effects under naive, ITT, CACE model ----

naive_itt_cace_estimation <- function(dv, label, scale, iv_robust=TRUE) { 
  # naive model based on observed treatment (disregarding assignment)
  observed_formula <- paste(dv, " ~ ", paste(c(covar_names, treatment), collapse=" + "))
  model_observed <- lm(observed_formula,  data = waves_df_wide)
  summary(model_observed)
  
  # ITT model based on encouragement assignment
  itt_formula <- paste(dv, " ~ ", paste(c(covar_names, assignment), collapse=" + "))
  model_itt <- lm(itt_formula,  data = waves_df_wide)
  summary(model_itt)
  
  # IV model to identify CACE
  iv_formula <- paste(dv, " ~ ", 
                      paste(c(covar_names, treatment), collapse=" + "), "|", 
                      paste(c(covar_names, assignment), collapse=" + "), sep = "") %>% as.formula
  if (iv_robust==FALSE){ model_cace <- ivreg(iv_formula,  data = waves_df_wide) }
  if (iv_robust==TRUE){ model_cace <- iv_robust(iv_formula,  data = waves_df_wide) }
  summary(model_cace, diagnostics = TRUE)
  
  # get estimates
  model_observed_predict <- tidy(model_observed)
  model_itt_predict <- tidy(model_itt)
  model_cace_predict <- tidy(model_cace)
  
  model_effects <- rbind(filter(model_observed_predict, str_detect(term, "wahlomat")),
                         filter(model_itt_predict, str_detect(term, "wahlomat")),
                         filter(model_cace_predict, str_detect(term, "wahlomat"))[
                           names(model_itt_predict)
                         ])
  
  model_effects$model <- c("Observed", "ITT", "CACE")
  model_effects$dv <- label
  model_effects$scale <- scale
  
  model_effects <- mutate(model_effects, 
                          lo95 = estimate - 1.96*std.error,
                          hi95 = estimate + 1.96*std.error,
                          lo83 = estimate - 1.4*std.error,
                          hi83 = estimate + 1.4*std.error)
  model_effects
}


### Effect heterogeneity ------------------------------------------

## Generic function to identify subgroup effects under naive, ITT, CACE model ----

subgroup_naive_itt_cace_estimation <- function(dv, label, scale, subgroup_var, labels_subgroup_var, iv_robust = TRUE, model_out = TRUE, estimand = "CACE") {
  
  # naive model based on observed treatment (disregarding assignment)
  observed_formula <- paste(dv, " ~ ", paste(c(covar_names, subgroup_varXtreatment), collapse=" + "))
  model_observed <- lm(observed_formula,  data = waves_df_wide)

  # ITT model based on encouragement assignment
  itt_formula <- paste(dv, " ~ ", paste(c(covar_names, subgroup_varXencouragement), collapse=" + "))
  model_itt <- lm(itt_formula,  data = waves_df_wide)

  # IV model to identify CACE
  iv_formula <- paste(dv, " ~ ", 
                      paste(covar_names, collapse = " + "), " + ",
                      paste(subgroup_varXtreatment, collapse=" + "), "|", 
                      paste(covar_names, collapse = " + "), " + ",
                      paste(subgroup_varXencouragement, collapse=" + "), sep = "") %>% as.formula
  if (iv_robust==FALSE){ model_cace <- ivreg(iv_formula,  data = waves_df_wide) }
  if (iv_robust==TRUE){ model_cace <- iv_robust(iv_formula,  data = waves_df_wide) }
  
  # get estimates
  model_observed_predict <- tidy(model_observed)
  model_itt_predict <- tidy(model_itt)
  model_cace_predict <- tidy(model_cace)
  
  model_effects <- rbind(filter(model_observed_predict, str_detect(term, "wahlomat")),
                         filter(model_itt_predict, str_detect(term, "wahlomat")),
                         filter(model_cace_predict, str_detect(term, "wahlomat"))[
                           names(model_itt_predict)
                           ])
  
  # compute effects
  if(estimand == "naive"){
    effects_subgroup <- c(
      model_observed_predict$estimate[which(model_observed_predict$term == "wahlomat_use_combined")],
      model_observed_predict$estimate[which(model_observed_predict$term == "wahlomat_use_combined")] + 
        model_observed_predict$estimate[(which(model_observed_predict$term == "wahlomat_use_combined")+1):length(model_observed_predict$term)])
    
    # compute effect se's
    effects_subgroup_ses <- c(model_observed_predict$std.error[which(model_observed_predict$term == "wahlomat_use_combined")], 
                              sqrt(diag(vcov(model_observed))["wahlomat_use_combined"] + 
                                     diag(vcov(model_observed))[(which(model_observed_predict$term == "wahlomat_use_combined") + 1):length(model_observed_predict$term)] + 
                                     2*vcov(model_observed)[which(model_observed_predict$term == "wahlomat_use_combined"), (which(model_observed_predict$term == "wahlomat_use_combined") + 1):length(model_observed_predict$term)]))
  }

  if(estimand == "ITT"){
    effects_subgroup <- c(
      model_itt_predict$estimate[which(model_itt_predict$term == "wahlomat_encouragement")],
      model_itt_predict$estimate[which(model_itt_predict$term == "wahlomat_encouragement")] + 
        model_itt_predict$estimate[(which(model_itt_predict$term == "wahlomat_encouragement")+1):length(model_itt_predict$term)])
    
    # compute effect se's
    effects_subgroup_ses <- c(model_itt_predict$std.error[which(model_itt_predict$term == "wahlomat_encouragement")], 
                              sqrt(diag(vcov(model_itt))["wahlomat_encouragement"] + 
                                     diag(vcov(model_itt))[(which(model_itt_predict$term == "wahlomat_encouragement") + 1):length(model_itt_predict$term)] + 
                                     2*vcov(model_itt)[which(model_itt_predict$term == "wahlomat_encouragement"), (which(model_itt_predict$term == "wahlomat_encouragement") + 1):length(model_itt_predict$term)]))  
  }
  
  if(estimand == "CACE"){
  effects_subgroup <- c(
    model_cace_predict$estimate[which(model_cace_predict$term == "wahlomat_use_combined")],
    model_cace_predict$estimate[which(model_cace_predict$term == "wahlomat_use_combined")] + 
      model_cace_predict$estimate[(which(model_cace_predict$term == "wahlomat_use_combined")+1):length(model_cace_predict$term)])
  
  # compute effect se's
  effects_subgroup_ses <- c(model_cace_predict$std.error[which(model_cace_predict$term == "wahlomat_use_combined")], 
                            sqrt(diag(vcov(model_cace))["wahlomat_use_combined"] + 
                                   diag(vcov(model_cace))[(which(model_cace_predict$term == "wahlomat_use_combined") + 1):length(model_cace_predict$term)] + 
                                   2*vcov(model_cace)[which(model_cace_predict$term == "wahlomat_use_combined"), (which(model_cace_predict$term == "wahlomat_use_combined") + 1):length(model_cace_predict$term)]))
  }
  
  # generate effect data frame
  model_effects <- data.frame(estimate = effects_subgroup,
                              std.error = effects_subgroup_ses)  
  model_effects$heterogeneity <- labels_subgroup_var
  model_effects$dv <- label
  model_effects$scale <- scale
  model_effects$subgroup_var <- subgroup_var
  model_effects <- mutate(model_effects, 
                          lo95 = estimate - 1.96*std.error,
                          hi95 = estimate + 1.96*std.error,
                          lo83 = estimate - 1.4*std.error,
                          hi83 = estimate + 1.4*std.error)
  if(model_out == TRUE){
  model_effects
  }else{
    if(estimand == "naive") {print(summary(model_observed))}
    if(estimand == "ITT") {print(summary(model_itt))}
    if(estimand == "CACE") {print(summary(model_cace, diagnostics = TRUE))}
  }
}



### By political interest -----------------------------

## Set up instrument, treatment, covariates
waves_df_wide$polinterest_1 <- car::recode(waves_df_wide$polinterest.1, "1:2=1;3=2;4:5=3")
waves_df_wide$polinterest_1[waves_df_wide$polinterest_1 == 6] <- NA
val_labels(waves_df_wide$polinterest_1) <- c(`not at all/\nrather not`= 1,
                                      `somewhat` = 2,
                                      `rather/\nvery` = 3)
assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "age.1", "educ.1", "hhincome.1", "leftright.1")
subgroup_var <- c("as_factor(polinterest_1)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("not at all/\nrather not", "somewhat", "rather/\nvery")

## Vote (reported
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot
model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)

dat <- filter(model_effects, heterogeneity == "somewhat")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("not at all/\nrather not", "somewhat", "rather/\nvery")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_political_interest.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)



### By civic knowledge -----------------------------

## Set up instrument, treatment, covariates
waves_df_wide$civic_knowledge_2 <- waves_df_wide$civic_knowledge.2
val_labels(waves_df_wide$civic_knowledge_2) <- c(`low`= 1,
                                                 `medium` = 2,
                                                 `high` = 3)

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "age.1", "leftright.1", "hhincome.1") # leave out pol interest and hhincome because of post-treatment bias
subgroup_var <- c("as_factor(civic_knowledge_2)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("low", "medium", "high")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot
model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)

dat <- filter(model_effects, heterogeneity == "medium")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("low", "medium", "high")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_civic_knowledge.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)



### By candidate recognition ---------------------

## Set up instrument, treatment, covariates
waves_df_wide$candidate_knowledge_2 <- car::recode(waves_df_wide$candidate_knowledge.2, "0:2=1;3:6=2;7:9=3")
val_labels(waves_df_wide$candidate_knowledge_2) <- c(`low`= 1,
                                                 `medium` = 2,
                                                 `high` = 3)

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "age.1", "leftright.1", "hhincome.1") # leave out pol interest and hhincome because of post-treatment bias
subgroup_var <- c("as_factor(candidate_knowledge_2)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("low", "medium", "high")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot
model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)

dat <- filter(model_effects, heterogeneity == "medium")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("low", "medium", "high")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_candidate_recognition.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)






### By education ---------

## Set up instrument, treatment, covariates

waves_df_wide <- filter(waves_df_wide, educ_1 < 4) # drop 24 "still in school" people

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "age.1", "leftright.1") # leave out pol interest and hhincome because of post-treatment bias
subgroup_var <- c("as_factor(educ_1)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("pre-sec'\nschool", "secondary\nschool", "high\nschool")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot
model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)

dat <- filter(model_effects, heterogeneity == "high\nschool")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity, levels = c("pre-sec'\nschool", "secondary\nschool", "high\nschool")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_education.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)



### By age groups -----------------------------

## Set up instrument, treatment, covariates
waves_df_wide$age_1 <- car::recode(waves_df_wide$age.1, "18:29='18-29'; 30:45='30-45'; 46:60='46-60'; 61:100='61+'")

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "leftright.1")
subgroup_var <- c("as_factor(age_1)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("18-29", "30-45", "46-60", "61+")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot
model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)

dat <- filter(model_effects, heterogeneity == "61+")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("18-29", "30-45", "46-60", "61+")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_age.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)



### By ideology -----------------------------

## Set up instrument, treatment, covariates

waves_df_wide$leftright_1 <- car::recode(as.numeric(waves_df_wide$leftright.1), "1:4='left-\nwing'; 5:7='center'; 8:11='right-\nwing'")

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("gender.1", "age.1", "educ.1", "hhincome.1")
subgroup_var <- c("as_factor(leftright_1)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("left-\nwing", "center", "right-\nwing")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot

model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)


dat <- filter(model_effects, heterogeneity == "center")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("left-\nwing", "center", "right-\nwing")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_ideology.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)




### By gender -----------------------------

## Set up instrument, treatment, covariates

waves_df_wide$gender_1 <- ifelse(waves_df_wide$gender.1 == 1, "male", "female")

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("age.1", "educ.1", "hhincome.1")
subgroup_var <- c("as_factor(gender_1)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("female", "male")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot

model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)


dat <- filter(model_effects, heterogeneity == "female")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("female", "male")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))
p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_gender.png",
       width = 3.5, height = 1.25, units = "cm", dpi = 300, scale = 5)




### By party preference -----------------------------

## Set up instrument, treatment, covariates

assignment <- "wahlomat_encouragement"
treatment <- "wahlomat_use_combined"
covar_names <- c("age.1", "educ.1", "hhincome.1")
subgroup_var <- c("as_factor(secondvote_3)")
subgroup_varXencouragement <- paste(subgroup_var, "wahlomat_encouragement", sep="*")
subgroup_varXtreatment <- paste(subgroup_var, "wahlomat_use_combined", sep="*")
labels_subgroup_var <- c("CDU/CSU", "SPD", "Left", "Greens", "FDP", "AfD", "Other", "DK")

## Vote (reported)
dv <- "voted_2017_yesno"
label <- "Turnout"
scale <- "[0-1]"
model_effects_voted_2017_yesno <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Change in vote choice
dv <- "secondvote_3to4"
label <- "Switched vote choice"
scale <- "[0-1]"
model_effects_secondvote_3to4 <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Knowledge on party positions
dv <- "wahlomat_perceivedpos_correct_avg"
label <- "Knowledge on parties'\nissue positions"
scale <- "[0-1]"
model_effects_wahlomat_perceivedpos_correct_avg <- subgroup_naive_itt_cace_estimation(dv, label, scale, subgroup_var, labels_subgroup_var)

## Combine model results; plot

model_effects <- rbind(model_effects_voted_2017_yesno,
                       model_effects_secondvote_3to4,
                       model_effects_wahlomat_perceivedpos_correct_avg)

dat <- filter(model_effects, heterogeneity == "SPD")
labels = paste(as.character(dat$dv), dat$scale)

model_effects <- transform(model_effects, 
                           heterogeneity = factor(heterogeneity,levels = c("CDU/CSU", "SPD", "Left", "Greens", "FDP", "AfD", "Other", "DK")),
                           dv = factor(dv, levels = unique(model_effects$dv), labels = labels))

p <- ggplot(model_effects, aes(heterogeneity, estimate)) + 
  geom_pointrange(aes(ymin = lo95, ymax = hi95), shape = 19, fatten = 1, size = 1) + 
  geom_pointrange(aes(ymin = lo83, ymax = hi83), shape = 19, fatten = 1, size = 2) + 
  geom_hline(yintercept = 0) + 
  facet_wrap( ~ dv, scales = "free", nrow = 1) +
  theme_ipsum_rc() + 
  theme(panel.background = element_blank(), 
        axis.title.y = element_blank(), 
        axis.title.x = element_blank(), 
        axis.text.x = element_text(size = 10, angle = 90), 
        plot.caption = element_text(size = 8), 
        strip.text.x = element_text(size = 12, colour = "darkblue", face = "bold"),
        plot.margin=unit(c(0.1,0.1,0.1,0.1), "cm"),
        panel.spacing = unit(0.25, "cm"))

p
# save plot
ggsave(p, file = "../output/wahlomat_effects_by_preference.png",
       width = 3.5, height = 1.75, units = "cm", dpi = 300, scale = 5)








